import numpy as np
import geopandas as gpd
import pandas as pd
import math
from shapely.geometry import MultiLineString,LineString,Polygon, Point,MultiPoint
import rtree

#region 生成匹配关系 
"""  找到gdf2中满足条件的line1，即gdf2：line1a，line1b,line1c....  """
def gen_matchrelation(gdf1, gdf2, buffer_size,arr,overlap_threshold=0.9):

    ##构建缓冲区
    buffer1 = gdf1.copy()
    buffer1['geometry'] = buffer1['geometry'].buffer(buffer_size, cap_style=2)
    buffer2 = gdf2.copy()
    buffer2['geometry'] = buffer2['geometry'].buffer(buffer_size, cap_style=2)
    # 构建树
    buffer1_idx = build_rtree_index(buffer1)
    buffer2 = buffer2.sort_index()
    match_list = []
    #遍历寻找
    for row in buffer2.itertuples():
        overlaps = buffer1_idx.intersection(row.geometry.bounds)
        lines_list = []
        for idx in overlaps:
            feature1 = buffer1.loc[idx]
            if row.geometry.intersects(feature1.geometry):
                line1 = gdf1.loc[gdf1['myid'] == feature1.myid].iloc[0]
                line2 = gdf2.loc[gdf2['myid'] == row.myid].iloc[0]
                overlap_lenth1 = line1.geometry.intersection(row.geometry).length
                overlap_lenth2 = line2.geometry.intersection(feature1.geometry).length

                # 匹配条件
                if (overlap_lenth1 / line1.geometry.length) > overlap_threshold or (overlap_lenth2 / line2.geometry.length) > overlap_threshold:
                    if  attribute_comparison(line1, line2,arr)==1:
                        lines_list.append(feature1.myid)

        if lines_list:
            lines_list.insert(0, row.myid)
        if lines_list:
            match_list.append(lines_list)

    df_lst = [pd.DataFrame([row]) for row in match_list]
    match_df = pd.concat(df_lst, axis=0, ignore_index=True)
    return match_df

def gen_matchrelationUA(gdf1, gdf2, buffer_size,overlap_threshold=0.9):

    ##构建缓冲区
    buffer1 = gdf1.copy()
    buffer1['geometry'] = buffer1['geometry'].buffer(buffer_size, cap_style=2)
    buffer2 = gdf2.copy()
    buffer2['geometry'] = buffer2['geometry'].buffer(buffer_size, cap_style=2)
    # 构建树
    buffer1_idx = build_rtree_index(buffer1)
    buffer2 = buffer2.sort_index()
    match_list = []
    #遍历寻找
    for row in buffer2.itertuples():
        overlaps = buffer1_idx.intersection(row.geometry.bounds)
        lines_list = []
        for idx in overlaps:
            feature1 = buffer1.loc[idx]
            if row.geometry.intersects(feature1.geometry):
                line1 = gdf1.loc[gdf1['myid'] == feature1.myid].iloc[0]
                line2 = gdf2.loc[gdf2['myid'] == row.myid].iloc[0]
                overlap_lenth1 = line1.geometry.intersection(row.geometry).length
                overlap_lenth2 = line2.geometry.intersection(feature1.geometry).length

                # 匹配条件
                if (overlap_lenth1 / line1.geometry.length) > overlap_threshold or (overlap_lenth2 / line2.geometry.length) > overlap_threshold:
                    lines_list.append(feature1.myid)

        if lines_list:
            lines_list.insert(0, row.myid)
        if lines_list:
            match_list.append(lines_list)

    df_lst = [pd.DataFrame([row]) for row in match_list]
    match_df = pd.concat(df_lst, axis=0, ignore_index=True)
    return match_df

""" 构建空间查询树 """
def build_rtree_index(gdf):
    gdf_idx = rtree.index.Index()
    for idx, geom in gdf.iterrows():
        gdf_idx.insert(idx, geom.geometry.bounds)
    return gdf_idx

""" 属性对比 """
def attribute_comparison(line1, line2,arr):
    a=0
    if line1[arr]==line2[arr]:
        a=1
    return(a)

""" 处理匹配关系的df """
def process_df(match_df):
    # 将NaN值替换为空字符串
    match_df.fillna('', inplace=True)  
    # 转置，去除重复行，再次转置
    match_df = match_df.T.drop_duplicates().T

    new_match_df = pd.DataFrame(columns=["col1", "col2"])
    # 遍历每一行
    for i, row in match_df.iterrows():
        col1 = row[0]
        for j in range(1, len(row)):
            # 获取第j列数据
            colj = row[j]
            if colj:  # 如果colj不为空
                new_row = {"col1": col1, "col2": colj}
                # 将新行添加到new_match_df中
                new_match_df = pd.concat([new_match_df, pd.DataFrame([new_row])], ignore_index=True)

    return new_match_df

###################################################################################################

#region  生成匹配要素之间的连接箭头 

def create_indictline_shapefile(line1,line2,df_matchresults,crs):
    lines = []
    ids1 = []
    ids2 = []
    a=0
    for row in df_matchresults.itertuples(): 

        for line_a in line1.itertuples():

            if line_a.myid == row.col1:
                
                a=a+1
                midpoint1 = line_a.geometry.centroid                
                Ppoint1 = line_a.geometry.interpolate(line_a.geometry.project(midpoint1))
                line_b=line2[line2['myid'] == row.col2]     
                if line_b.empty:
                    continue
                midpoint2 = line_b.geometry.centroid
                Ppoint2 = line_b.interpolate(line_b.project(midpoint2))

                line=LineString([(Ppoint1.x,Ppoint1.y), (Ppoint2.x,Ppoint2.y)])
                lines.append(line)
                ids1.append(row.col1)
                ids2.append(row.col2)

    print(a)
    d = {'col1': ids1, 'col2': ids2}
    ID = pd.DataFrame(d)
    new_gdf = gpd.GeoDataFrame(ID, geometry=lines, crs=crs)
    return new_gdf

###################################################################################################
#region  给csv增加shp的属性 

def add_shplen(line1,line2,df):
    """ 给csv增加shp的长度比 """
    line_id='myid'
    
    # lengths1 = df.apply(lambda row: 
    #     (print(f"当前循环的myid: {row[0]}"),  # 打印当前myid
    #      line1.loc[line1[line_id] == str(row[0]), 'geometry'].iloc[0].length)[1], 
    #     axis=1)

    lengths2 = df.apply(lambda row: 
        (print(f"当前循环的myid: {row[-1]}"),  # 打印当前myid
         line2.loc[line2[line_id] == str(row[-1]), 'geometry'].iloc[0].length)[1], 
        axis=1)   
    # Lena=lengths1
    # lengths1 = df.apply(lambda row: line1.loc[line1[line_id] == str(row.col1), 'geometry'].iloc[0].length, axis=1)
    # lengths2 = df.apply(lambda row: line2.loc[line2[line_id] == str(row.col2), 'geometry'].iloc[0].length, axis=1)
    LR=lengths2

    # Concatenate the lengths Series with the original DataFrame
    df = pd.concat([df, LR.rename('LR')], axis=1)
    return df


def add_Overlap(line1, line2, df, buffer_radius=100):
    """ 给csv增加叠置面积比 """
    def calculate_overlap(row):
        # 提取当前行的ID
        id1 = row['col1']
        id2 = row['col2']
        
        # 获取对应几何对象
        geom1 = line1[line1['myid'] == id1].geometry.iloc[0]
        geom2 = line2[line2['myid'] == id2].geometry.iloc[0]
        
        # 构建缓冲区
        buffer1 = geom1.buffer(buffer_radius)
        buffer2 = geom2.buffer(buffer_radius)
        
        # 计算叠置面积
        intersection = buffer1.intersection(buffer2)
        overlap_area = intersection.area if not intersection.is_empty else 0
        
        # 计算最大面积占比
        area1 = buffer1.area
        area2 = buffer2.area
        min_area = min(area1, area2)
        ratio = overlap_area / min_area if min_area > 0 else 0
        
        return pd.Series([overlap_area, ratio], index=['Overlap', 'Overlap_Ratio'])

    # 应用计算并合并结果
    df = df.merge(
        df.apply(calculate_overlap, axis=1),
        left_index=True,
        right_index=True
    )
    return df

###################################################################################################


#region 将等高线划分为开放式/闭合式等 
def classify_contours(contour1):

    boundary_contour1 = contour1.geometry.unary_union.envelope.boundary
    contour1['gcla'] = 3
    for idx, row in contour1.iterrows():
        geom = row.geometry
        ## 识别闭合线
        if geom.is_closed:
            contour1.at[idx, 'gcla'] = 2
            continue

        # 非闭合情况下计算起点终点距离
        coords = list(geom.coords)
        start = Point(coords[0])
        end = Point(coords[-1])
        distance_start_end = start.distance(end)

        # 计算最小旋转外接矩形及其短边长度
        mrr = geom.minimum_rotated_rectangle
        if isinstance(mrr, Polygon):
            mrr_coords = list(mrr.exterior.coords)
            edge0 = Point(mrr_coords[0]).distance(Point(mrr_coords[1]))
            edge1 = Point(mrr_coords[1]).distance(Point(mrr_coords[2]))
            short_side = min(edge0, edge1)/10
        else:
            short_side = 0

        if distance_start_end < short_side:
            contour1.at[idx, 'gcla'] = 2
            continue

        ##识别开放线     
        buffer_area1 = boundary_contour1.buffer(100) 
        buffer_area2 = boundary_contour1.buffer(5) 

        # intersection = geom.intersection(buffer_area)
        # intersectionP = geom.intersection(boundary_contour1.buffer(5) )


        if not geom.within(buffer_area1) and start.within(buffer_area2) and end.within(buffer_area2) :
            contour1.at[idx, 'gcla'] = 1





    # 保存结果
    return contour1


###################################################################################################






























""" 计算弯曲度 """
def Cal_ST(geodf):

    # 遍历GeoDataFrame，计算每个LineString的弯曲度
    for idx, row in geodf.iterrows():
        geo = row.geometry
        length=row.geometry.length
        if geo.is_closed:
            L= geo.length
            coords = list(geo.coords)
            polygon = Polygon(coords + [coords[0]])
            # 计算多边形的长边长度
            mrrect = polygon.minimum_rotated_rectangle
            corners = mrrect.exterior.coords  # 获取矩形四个顶点的坐标
            x_diff =  Point(corners[1]).distance(Point(corners[2]))
            y_diff =  Point(corners[1]).distance(Point(corners[0]))
            st=L/(2*(x_diff+y_diff))


        else:
            start_point = geo.coords[0]
            end_point = geo.coords[-1]
            dis = LineString([start_point, end_point]).length
            L= geo.length 
            st=dis/L
            # 计算弯曲度并保存到相应的列中
        geodf.at[idx, 'R1ST'] = st
        geodf.at[idx, 'L'] = length
   
    return geodf

################################################################################
""" 计算邻接要素数量 """
def CalAdjacentNum(geodf, distance=0.1):   

    sindex = geodf.sindex
    geodf['adn'] = 0

    # 遍历每个要素，计算邻接数量
    for idx, row in geodf.iterrows():
        adn_count = count_adjacent(row.geometry, geodf,sindex, distance)
        # 保存结果
        geodf.at[idx, 'adn'] = adn_count

    return geodf

def count_adjacent(geom, geodf,index, distance=0.1):

    buffer_geom = geom.buffer(distance)
    intersecting_ids = list(index.intersection(buffer_geom.bounds))
    # 检查每个相交的要素是否真的与缓冲区相交
    adjacent_count = 0
    for idx in intersecting_ids:
        adjacent_geom = geodf.geometry[idx]
        if buffer_geom.intersects(adjacent_geom):
            adjacent_count += 1
    adjacent_count=adjacent_count-1        
    return adjacent_count